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Individual-based models for bacterial chemotaxis in the diffusion 

aymptotics 

Mathias Rousset* and Giovanni Samaeyt 

Wc discuss velocity-jump models for chemotaxis of bacteria with an internal state that 
allows the velocity jump rate to depend on the memory of the chemoattractant concen- 
tration along their path of motion. Using probabilistic techniques, we provide a pathwisc 
result that shows that the considered process converges to an advection-diffusion pro- 
cess in the (long-time) diffusion limit. We also (re-)prove using the same approach that 
the same limiting equation arises for a related, simpler process with direct sensing of 
the chemoattractant gradient. Additionally, we propose a time discretization technique 
that retains these diffusion limits exactl y, i .e.. without error that depends on the time 
discretization. In the companion paper 121J these results are used to construct a cou- 
pling technique that allows numerical simulation of the process with internal state with 
asymptotic variance reduction, in the sense that the variance vanishes in the diffusion 
limit. 
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1. Introduction 

Generally, the motion of flagellated bacteria consists of a sequence of run phases, 
during which a bacterium moves in a straight line at constant speed. The bac- 
terium changes direction in a tumble phase, which takes much less time than the 
run phase and acts as a reorientation. To bias movement towards regions with high 
concentration of chemoattractant, the bacterium adjusts its turning rate to increase, 
resp. decrease, the prob ability of tumbling when moving in an unfavorable, resp. fa- 

lOIOQI 

vorable, direction 1 1 I Since many species are unable to sense chemoattractant 
gradients reliably due to their small size, this adjustment is often done via an intra- 
cellular mechanism that allows the bacterium to retain information on the history 
of the chemoattractant concentrations along its path SI. The resulting model can be 
formulated as a velocity-jump process, combined with an ordinary differential equa- 
tion (ODE) that describes the evolution of an internal state that incorporates this 
memory effect UJ2i . This model will be called the "fine-scale" or "internal dynamics" 
model in this paper. 

The probability density of the velocity-jump process evolves according to a ki- 
netic equation, in which, besides position and velocity, the internal variables appear 
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2 Individual-based model of chemotaxis 

as additional dimensions. A direct deterministic simulation of the density distribu- 
tion of all the variables of the model is therefore prohibitively expensive. Hence, 
it is of interest to study the relation of this model with simplified, coarse-level 
descriptions of the bacteria dynamics. 

In chemotaxis, a first coarse description is obtained by neglecting velocity and 
internal variables and considering the bacterial position density on large space and 
time scales. One then postulates an advection-diffusion equation for this bacterial 
position density, in which a chemotactic sensitivity coefficient incorporates the effect 
of chemoattractant concentrations gradients on the dens ity flux es. This assumption 
leads to the classical Keller-Segel equations (see and I 13 | 14 | f or numerous histor- 
ical references) . Several works have also considered the motion of a bacterium to be 
governed by a veloc ity- jump process corresponding to a kinetic description for the 
phase-space density HOLlEST j n these models, the velocity jump rate of the bacteria 
depends on the local chemoattractant gradient. In the present paper, the model 
corresponding to such bacteria will be called a "coarse" or "d irect gr adient sens- 
ing" model. Unlike the limiting Keller-Segel equation, see e.g. I 5 | H | 1^ I ; t ne kinetic 
description does not exhibit finite-time blow-up (which is believed to be unphysi- 
cal) under certain biologically relevant assumptions on the turning kernel^. These 
models also converge to a Keller-Segel equation in the appropriate drift-diffusion 
limit, see e.g. KHUHS InlSl, the kinetic equation associated to a model with internal 
state (similar to the model considered here) has also be shown to formally converge 
to a Keller-Segel equation using moment closure assumptions and an appropriate 
(diffusion) scaling; a similar result using an infinite moment expansion has been 
obtained in I^SI. Then, the para mete rs of internal dynamics appear in the expression 
for the chemotactic sensitivity EE Existence results when the model is coupled to 
a mean-field production of chemo-attractants are also available ^, and analysis of 
the long-time behavior is carried out in Moreover, some recent studies have 
assessed the biological relevance of such a model through experimental validation 
of travelling pulses ^2. 

The present paper contributes to the analysis of velocity-jump processes for 
chemotaxis of bacteria with internal dynamics, and their numerical discretization, 
in the diffusive asymptotics. The main contributions are twofold : 

• From an analysis point of view, we consider (fine-scale) velocity-jump pro- 
cesses with internal dynamics (|2.5|) . with jump rates satisfying (|2.4j) and 
internal dynamics satisfying the two main assumptions (|3.2[) and (|3.3[) as 
well as other technical assumptions contained in Section 13.31 We rigor- 
ously prove, using probabilistic arguments, the pathwise convergence of the 
velocity-jump process described above to the stochastic differential equa- 
tion (SDE) (|4.1I) in the long time limit (see Proposition ^. lip . The technical 
steps of the proof are given in the beginning of Section We also recall a 
convergence result for the (coarse) process with direct gradient sensing, see 
Proposition 14.31 For an appropriate choice of the model parameters, both 
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limits can be made identical. 
• From a numerical point of view, we introduce a time-discretized model (see 
Section [5]) for both the fine-scale and the coarse velocity-jump processes. 
When discretizing the fine-scale process with internal dynamics, particular 
attention is paid to the time discretization errors due to the approximate 
solution of the ODE describing the evolution of the internal state. As a 
result, the asymptotic analysis of the time-continuous case still holds in 
an exact fashion, ensuring that the advection-diffusion limits are preserved 
after time discretization. 

Due to the possibly high number of dimensions of the kinetic model with inter- 
nal state, the evolution of the bacterial density away from the diffusive limit may 
rather be simulated using a stochastic particle method. However, a direct stochastic 
particle-based simulation suffers from a large statistical variance, even in the diffu- 
sive asymptotic regime, raisi ng the important issue of variance reduced simulation. 
In the companion paper ^H, we will use the results of the analysis presented here 
to construct a numerical scheme that couples the fine-scale model for bacteria with 
internal dynamics and the simpler, coarse model for bacteria with direct gradient 
sensing to obtain a variance-reduced simulation of the fine-scale process. We show 
that the variance reduction is asymptotic, in the sense that the statistical variance 
scales as the error between the fine-scale and coarse description. 

This paper is organized as follows. In Section [2j we discuss the velocity-jump 
processes that will be considered, introducing both the fine-scale process with inter- 
nal state and a simpler, coarse process in which the internal state has been replaced 
by direct gradient sensing. We discuss the relevant asymptotic regimes in Section[3l 
and analyze the limits of the considered processes in these regimes in Section 0] 
Section [5] discusses a discretization scheme for the velocity-jump processes that re- 
tains the diffusion limit of the continuum processes. Section [6] contains conclusions 
and a more detailed discussion of the relation with the companion paper ^H. 



2. Particle-based models for bacterial chemotaxis 
2.1. Model with internal state 

We consider bacteria that are sensitive to the concentration of m chemoattractants 
{pi{x))^L 1 , with pi(x) >0 for x£R d . While we do not consider time dependence of 
chemoattractants via production or consumption by the bacteria, a generalization 
to this situation is straightforward. Bacteria move with a constant speed v (run), 
and change direction at random instances in time (tumble), in an attempt to move 
towards regions with high chemoattractant concentrations. As in we describe this 
behavior by a velocity-jump process driven by some internal state y € Y C K" of each 
individual bacterium. The internal state models the memory of the bacterium and 
is subject to an evolution mechanism attracted by a function ijj : K m — > M. n of the 
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chemoattractants concentrations, 

S(x):=ii Pl (x),...,p m (x))eR n , 

where x is the present position of the bacterium. A typical choice is 2m = n, 
ip{pi,. . .,p m ) — (pi, - ■ ■ ,Prn,0, ■ ■ ■ ,0), and the internal mechanism being such that the 
TTi first internal variables are memorizing the concentrations of chemoattractants, 
while the m last internal variables are computing the variations of concentrations 
along a bacterium trajectory (see Section |2~2|) . S and its two first derivatives are 
assumed to be smooth and bounded. 

We use (l s ,t\) as a reference for length and time, where 

• l s is the typical length of the chemoattractant concentration variations 
(which we assume similar for all the different species of chemoattractants). 

• t\ is the typical time between two changes of the bacterium velocity direc- 
tion (tumbling). 

We also introduce the time t x that indicates the typical time on which a bacterium 
observes changes in chemoattractant concentration; the typical speed of a bacterium 
is thus l s /t x . Anticipating the time-scale analysis in Section[3l we already introduce 
a small parameter 

{eq:eps} e:=^«l, (2.1) 

indicating that we assume the typical time between two velocity changes to be much 
smaller than the typical time on which we observe macroscopic movement of the 
bacteria. 

The model is now presented directly in dimensionless form. The evolution of 
each bacterium individual position is denoted as 

with velocity 

at 

where S d_1 represents the unit sphere in R d . The evolution of the internal state is 
denoted by t^Y t . The internal state adapts to the local chemoattractant concen- 
tration via an ordinary differential equation (ODE), 

{eq:internal} -±=F e (Y t ,S(X t )), (2.2) 

at 

which is required to have a unique fixed point y* = S(x*) for every fixed value x* S 
K d . We also introduce the deviations from equilibrium z — S(x)—y. The evolution 
of these deviations is denoted as 

t^Z t = S(X t )-Y t . 
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The velocity of each bacterium is switched at random jump times (T ra ) n >i that are 
generated via a Poisson process with a time dependent rate given by X(Z t ), where 
z i — y A(z) is a smooth function satisfying 



ratebound} < A m ; n < A (z) < A max . 

We will denote its expansion for small z as 

\(z) = \ Q -b T z + c x o(\z\ k ) 

with fc>2, &eR ra ; the constant c\ is used to keep track of the nonlinearity of A in 
the analysis. 

The new velocity V n at time T n is generated at random according to a centered 
probability distribution A4(dv) with J v A4(dv) =0, typically 

A4(dv) = a S d-i (dv), 

where <7§d-i is the uniform distribution on the unit sphere. 

The resulting fine-scale stochastic evolution of a bacterium is then described by 
a left continuous with right limits (lcrl) process 

t^(X t ,V t ,Y t ), 

which satisfies the following differential velocity-jump equation: 
(dX t 



(2.3) 

(2.4) {eq:lin_rate} 



dt 



-^=F e (Y t ,S(X t )) 



(2.5) {eq:process_noscale} 



x(z t )dt=e n+1 , with z t :=s(x t )-Y t 

Vt=V n forie[T„,T n+1 ] , 

with initial condition X ,V eR d , Y eW l and T o = 0. In ([23]), (6 n ) n>1 denote i.i.d. 
random variables with normalized exponential distribution, and (V n ) n>1 denote 
i.i.d. random variables with distribution M.(dv). 

Remark 2.1. In one spatial dimension, assuming that A4(dv) = \ (6 + 1 +<5_i), the 
system (|2.5[) is equivalent to the following system: 



dXt 
dt 

dY L _ 
dt 



T„ 



--eV t 

F e (X u S(Y t )) 
1 \{Z t )dt = 26 n+1 , 



(2.6) {eq:processld} 



V t = V n for t G [T n ,T n+1 ] with V n = -V n _ 



I 
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in the sense that both processes have the same probability path distribution. This 
can be checked on the Markov infinitesimal generator of probability transitions, 
see (|3.5|) and discussion below. (Note that, in this process, the velocity is always 
reversed. This is compensated by a factor 2 in the equation defining the jump times.) 



2.2. Example 

[71 

For concreteness, we provide a specific example, adapted from L H We consider m = 1, 
i.e. there is only one chemoattractant p\{x). We first describe a cartoon dynamics ex- 
hibiting an excitation-adaptation behaviour. The internal state is two-dimensional, 
i.e. n = 2, and y= (2/1,3/2) satisfies the following ODE: 

dj/ift) = pi(x)- yi (t) 

e: exc-ad} i * *° (2.7) 

^ dy 2 (t) _ -p 1 (x)+y 1 (t)-y 2 (t) 

dt t e 

in which t a , resp. i e , represent an adaptation, resp. excitation, time; and yi, 
resp. 2/2, represent the adapting, resp. exciting, variable. This model has a single 
fixed point y* = (y^,y^) = (p 1 (x),0) =: S(x), and the deviation variables are given by 
z=(z\,z%) = S{x) — y= (pi{x) — 2/2)' When the excitation time is much smaller 
than the adaptation time, t e <^it a) the variable y\ therefore adapts "slowly" to the 
environment and memorizes it, while 2/2 computes "faster" the lag y\ — pi(x), giving 
the response to changing environments. In the asymptotics t e /t a —>0, the internal 
dynamics reduce to a scalar equation 

scalar-yj — = , (2.8) 

at t a 

while the exciting variable reduces to the difference 2/2 = Vi — Pi{x) = —z\ which is 
instantly learned by the bacteria. In this limit, we may eliminate the second internal 
state variable, and write y = y\ and Six) :=p±(x), and z = S(x) — y. 

For the turning rate z<— >\(z), we choose a nonlinear strictly decreasing smooth 
function, depending on a scalar £ that aims at measuring the difference between 
the chemoattractant concentration and the memory variable. When considering the 
two-dimensional internal dynamics, we choose C = z 2 = — IJ2, while, in the scalar case 
C = z. As a consequence, C>0, resp. C<0, give the information to the bacterium 
that it is moving in a favorable, resp. unfavorable, direction. A specific choice of 
A(C) is 

{eqirate} A(C) = 2A Q - ^ arctan U^/CJ J , /3>0. (2.9) 



Note that this turning rate satisfies 

A(z)-A -/?C + o(|C| 3 ) 



I 
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2.3. Model with direct gradient sensing 

We now turn to a simplified model, in which the internal process (|2.2j) . and the 
corresponding state variables, are eliminated. Instead, the turning rate depends 
directly on the chemoattractant gradient. The process with direct gradient sensing 
is a Markov process in position- velocity variables 

t^(x?,v t % 

which evolve according to the differential velocity-jump equations 

■ eVf 



AX? 



dt 



K{xiv t c y\t=e n+l 



(2.10) {eq: cprocess} 



V t c = V n forie[T n c ,T n c 



+u 



denote i.i.d. random variables 



with initial condition X ,V eR d . In (j2~T0)) . 
with normalized exponential distribution, and (V„) Jl>1 denote i.i.d. random vari- 
ables with distribution M(dv). The turning rate of the process with direct gradient 
sensing is assumed to be of the form 



\ c e (x,v)~\ Q ~eAj(x)v + 0(e 2 ), 



and to satisfy 



< A min < X(x,v) < X 



max- 



(2.11) {eq: control_rate} 



(2.12) {eq:rateboundc} 



The vector field A e (x) S R d is usually a linear combination of the columns of VS(x) £ 
R dx ™, and, as a consequence on the smoothness assumption on S(x), satisfies that, 
for e — 0, both A e and its derivatives approach a limiting vector field Aq in uniform 
norm. 



The model (|2 . 1 1 1) may describe a large bacterium that is able to directly sense 
chemoattractant gradients. In the case n — 1, the turning rate (|2 . 1 1 1) is proportional 
to VS(x)v £ R, which can be interpreted as follows: the rate at which a bacterium 
will change its velocity direction depends on the alignment of the velocity with 
the gradient of the chemoattractant concentration \7S(x), resulting in a transport 
towards areas with higher chemoattractant concentrations. 



3. Asymptotic regimes and kinetic formulation 
3.1. Notation 

Throughout the text, the Landau symbol O denotes a deterministic and globally 
Lipschitz function (with O(0) = 0). Its precise value may vary from line to line, and 
may depend on the parameters of the model. However, its Lipschitz constant is 
uniform in e. In the same way, we will denote generically by C > a deterministic 
constant that may depend on all the parameters of the model except e. 
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3.2. Definition of time scales 

Recall that we have introduced two time scales, characterized by (i) the typical 
time t x on which a bacterium observes changes in chemoattractant concentration; 
and (ii) the typical time t\ between two changes of the bacterium velocity direction 
(tumbling). We have assumed a large time-scale separation between these two time 
scales by introducing the small parameter e = t\/t x ^l. 

Let us now turn to the typical time t v associated with the dynamics of the 
internal state. In principle, several typical time scales (ti, v 3>i2.j/ 3> •••) may be 
required to describe the evolution of the internal state; t y :=t\ >y is assumed to be 
the largest one. The vector field driving the ODE (12.21) of the internal state is thus 
of order 

tx 



'- li 



and we will assume that 



{eq:F} F e ~^ = e 1 - ff , (3.1) 



ty 

with 5>l/k, where k is defined by (|2.4|) . 

Consequently, we will consider three time scales : 

• A fast time scale of order 0(1) given by the rate of change of the velocity 
direction. 

• A time scale of order C(l/e 1_<5 ) with 5>l/k given by the internal state 
evolution; depending on S, this time scale can be very fast, fast or interme- 
diate. 

• A slow time scale 0(1/ e) given by the evolution of the chemoattractant 
concentration as seen from the bacteria. 

3.3. Asymptotic assumptions 

The loose assumption on the scale, equation (|3.1[) for the internal dynamics (|2.2|) . 
needs to be made precise by the following set of assumptions. First, we assume that 
the ODE (12. 2|) driving the internal state is well approximated by a near-equilibrium 
evolution equation : 

Assumption 3.1. We have 

{eq : tau_eps} F e (y, S ) = -T- 1 (y-s) + e 1 - 5 c F 0(\s-y\ 2 ), (3.2) 

where r e eR nx " is an invertible constant matrix. The constant cp is used to keep 
track of the dependence on the non-linearity of F t in the analysis. 

We now formulate the necessary assumptions on the scale of the matrix r e : 
Assumption 3.2. There is a constant C > such that for any t > 0, one has 



{eq: tau_scales} 



e -*r e 



<Ce-^ S / c . (3.3) 
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Assumption 13.21 is used to ensure exponential convergence of the linear ODE 
with time scale at least of order e 1- * 5 . When r e _1 is symmetric, Assumption 13.2 
follows from the assumption that the eigenvalues of r e _1 are strictly positive, with 
a lower bound of order 0(e 1 ^ 5 ). 

Assumption 3.3. There is a constant C > such that 



sup 

t>o 



<Ce~\ 



Assumption 13.31 is necessary for technical reasons in Lemma 14.81 When r e _1 is 
symmetric, then Assumption 13.31 follows from the assumption that the eigenvalues 
of t £ -1 are non-negative. 

Finally, we assume that the solution of the ODE driving the internal state in 
T5"|) satisfies the following long time behavior : 



Assumption 3.4. Consider a process in- St such that sup tg [ +oc ] \dSt/dt\ = 0(e), 
and assume \S -Y \ = 0(e s ). Then the solution of 



satisfies 



sup \Y t -S t \ = 0(e s ). 



te[o,+oo] 

Assumption 13.41 is motivated by the case of linear ODEs; indeed, in that case, 
Assumption 13.41 is a consequence of Assumption 13. 2\ as is seen from the following 
Lemma : 

Lemma 3.5. Assume that 

Fe(y,s) = -T e - 1 {y-s), 
then Assumption \3.2\ implies Assumption\3.4\ 



Proof. Denoting Z t ~ St — Y t in Assumption 13.41 Duhamel's integration formula 
yields 

so that by Assumption 13. 2\ we have 

\Zt\<\Z \+Ce f e-^'V^'c&t'. 



^iZol+Cexe^l-e-^^), 
and Assumption 13.41 follows. □ 



In a similar fashion, Assumption 13.31 implies the following technical bound 
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Lemma 3.6. Consider a process tt-tSt such that sup te j +oo] \dSt/dt\ = 0(e). As- 
sume that \Sq — Yq\ = 0(1), and that Assumvtion \3.l[ Assumption \37^ Assump- 
tion\S.3\ and Assumption\3.4\hold. Then, we have 



sup 

t€[0,+o 



(St-Yt) = 0(1). 



(3.4) {eq: techbound} 



Proof. Let us denote A t := r e 1 e T ' t (St — Yt), which satisfies 

dA t „ _1 , , -1 -r," 1 * Z' ^St 



Here again, Duhamel integration yields 



-F t (Y u S t )+T~\S t -Y t 



From Assumption EH] and Assumption [331 1-^(1*',-%') -T e -1 (St' -Y t >)\ — 
0(e 1+s )<0{e) so that 



\A t \<A - 



<A - 



0(e)di' 



-te y -*/C 



)0(e) 



<0(l) + 0(e- 1 )0(e) = 0(l), 
where we have used in the above Assumption 13.21 and Assumption 



□ 

In the remainder of the paper, we will assume that the solution of (|2.5p satisfies 

\Zo\ = 0(e s ), 

as well as: 

1^-^01=0(1), 
which implies, according to Assumption 13.41 and Lemma [ 



that 



sup {\Z t \) = 0{e% 
*e[o,+oo] 



as well as 



sup | 

te[o,+oo] 



= 0(1). 



Note that all the assumptions of the present section hold in the following case: (i) 
the ODE is linear (F t is linear), (ii) the associated matrix r e is symmetric positive 
and satisfies 



for some 8> 0, (hi) the initial internal state is close to equilibrium in the sense that 
t c ~ 1 Zq = 0(1). The assumptions of the present section can be seen as technical 
generalizations to non-linear ODEs with non-symmetric linear part. 
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3.4. Kinetic formulations 

The probability distribution density of the process with internal dynamics at time 
t with respect to the measure dxAi{dv)dy is denoted as p(x,v,y,t), suppressing 
the dependence on e for notational convenience, and evolves according to the Kolo- 
mogorov forward evolution equation (or master equation). In the present context, 
the latter is the following kinetic equation 

d t p + ev ■ V xP + div y (F e (y,S(x))p) = A (S(x) - y) (R(p) -p), (3.5) {e : kinetic} 

where 

R(p):= f p(;V,-)M(dv) 

is the operator integrating velocities with respect to M.. Similarly, the distribution 
density of the process with direct gradient sensing at time t with respect to the mea- 
sure dxA4(dv) is denoted as p c (x,v,t), and evolves according to the kinetic/master 
equation 

d t p c + ev- V x p c — (R(\ c t (x,v)p) — \ c e (x,v)p) . (3.6) {e :kinetic_control} 

We refer to ^ fo r the derivat ion of master equations associated to Markov jump 
processes, and to | 12 | 7 | 26 | 6|18 | f Qr an aS y m ptotic analysis using moment closures. 

4. Asymptotic analysis of individual-based models 

In this section, we investigate in more detail the asymptotic behavior of the pro- 
cess with internal dynamics and the simplified process with direct gradient sensing. 
We rigorously prove pathwise convergence to an advection-diffusion limit of both 
the process with direct gradient sensing (Section 14.11) and the process with inter- 
nal dynamics (Section I4.2|) . Choosing the model parameters appropriately, both 
advection-diffusion limits can be made identical. The two proofs are based on an 
asymptotic expansion of the jump times (the time between two tumble phases), and 
a comparison with a simpler random walk. They can be split into 4 steps. 

Step (i) Consider the evolution of the bacterium indexed by the number of tum- 
ble phases (jump times), as well as a proper Markovian random walk 
approximating the former. Show tightness of the coupled sequence (in 
fact their time continuous version obtained by linear interpolation). 
Step (ii) Show, using a Gronwall argument, that the two processes have the same 
advection-diffusion limit. This step relies on an asympotic expansion of 
the jump times in e. 

Step (iii) Compute the advection-diffusion limit of the Markovian random walk. 
Step (iv) Use a random time change to transfer the result to the bacteria position 
process itself. 
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4.1. Diffusion limit of the process with direct gradient sensing 

We first derive a pathwise advection-diffusion limit of the process with direct gradi- 
ent sensing (|2.10[) when e — >■ using classical probabilistic arguments. It should be 
noted that this diffusive limit as already been extensively studied in the context of 
bacterial chemotaxis, see, e.g. for justifications of Hilbert expansions at the PDE 
level, and ^1 where the diffusive limit of some more general class of stochastic pro- 
cesses is suggested. However, we give the probabilistic arguments in extenso since 
they will be of interest for the analysis of both the model with internal state, and of 
the co upl ing in the variance-reduced numerical scheme, discussed in the companion 
paper 1=^4 Moreover, they do not seem to appear as such in the literature. 
Diffusive times are denoted by 



{eq: cprocess_hydro} 



and the process with direct gradient sensing considered on diffusive time scales by 



\rc,e v c 
A t — A t/e 



This section shows that this process converges for e — > towards an advection- 
diffusion process, satisfying the stochastic differential equation (SDE) 



'DMxf) 

A 



1/2 



(4.1) 



{eq:maxw} 



{eq: cdensity} 



{eq: cdens_hydro} 



where 1 1-> is a standard Brownian motion, and the diffusion matrix is given by 
the covariance of the Maxwellian distribution: 



D = 



v®vM(dv)&. 



ndxd 



(4.2) 



In particular, this result implies at the level of the Kolomogorov/master evolution 
equation, that the position bacterial density 



n c - e (x,i):=n c (x,t/e 2 ):= / p c {x,v,t/e 2 )M(dv) 



satisfies the advection-diffusion equation 

1 



dm c >°: 



An 



div x (DV x n c '° - DA (x)n 



(4.3) 



(4.4) 



on diffusive time scales as e— >0. 

The proof relies on a perturbation analysis of the jump time differences, which 
is given by the following lemma : 

Lemma 4.1. The difference between two jump times of the process with direct 
gradient sensing satisfies 



{eq:deltaTc} 



•— 1 n+l 1 n 



Jn+1 

"a7 



+ e 



'n+l aT 



A§ 



Ai{X c Tc )V n + 9 n+l O{e 2 ). 



(4.5) 
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Proof. Using (|2.11l) . we get 

K{X c tl V t c ) = A - eA T t {X c t )V n + 0(e 2 ) 

= X -eAj(X^)V n + (\\\7A\\ oo + l)O(e 2 ) 

= X -eAj(X^)V n + O(e 2 ), 

where we have used differentiability of A e (x), and the property ||AT t c — Xj, c \\ <e(t- 
T£). Integrating on the time interval [T^,T^ +1 ] yields 

6 n+1 = AT c n+1 (Ao - eAj(X^ )V n + C3(e 2 )) , 
AT„ C +1 = - 



so that 



y n +i 



X -eAf(X^)V n + O(e 2 ) 

'n+1 / 1 



Ao 



1 



e— Ai (X^)V n )+9 n+1 0(e'). 

A 



□ 



We also need the following technical lemma : 

Lemma 4.2 (Random time change). Let (?H>a e (i)) e>0 a sequence of contin- 
uous, strictly increasing random time changes in R + , which, for e — > 0, converges 
in distribution (on any time interval for the uniform convergence topology) towards 
the identity function it- >i. Let (t—^X € (t)) e>0 a sequence of continuous random pro- 
cesses, such that 



t^-x: 



.(*) 



converges in distribution when e — > (on any time interval for the uniform 
convergence topology) towards a continuous process i—^X^. Then the sequence 
(t— >X^)^ >Q also converges in distribution when e— >0 (on any time interval for 
the uniform convergence topology) towards i—tX®. 



Proof. Since by assumption (£~i— >-a £ (F)) c >o converges to a deterministic limit, the 
coupled process tt-> (X^ a\,ct e (i}) also converges in distribution to (X®,t), and 
we can consider a Skorokhod embedding (a new probability space) associated with 
a e-sequence where this convergence holds almost surely. Then we can write (A 
denotes the minimum, and 1 the identity function): 



sup \X§-X9\M 
te[o,T] 



< sup 

te[o,T] 



XI 



A1+ sup 

te[o,T] 



< sup 

5£[0,a7 X (T)] 



x:- 



< sup 

se[0,2T] 



X 



<*•(«) 



x° 



A 1 + sup 

se[0,ar 1 (T)] 



A 1 + sup 

sS[0,2T] 



(*) 



Al, 



00' 



Al 



Al + 2 X 1 Q =(2T)<T 
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and the right hand side converges almost surely to using the (uniform) continuity 

Then the following result holds : 

Proposition 4.3. Consider the gradient sensing process described in Section \2.S\ 
The process i <— >X^' e = X^ e2 converges in distribution (with respect to the topology 
of uniform convergence) towards it-^X^°, solution of the SDE As a conse- 

quence, the density n c,e of the process with direct gradient sensing defined in (|4.3[) 
converges weakly towards n c, ° , solution of fr4-4\ )- 

Proof. The proof requires the definition of some auxiliary processes. First, we 
perform a random time change on the process of positions of the bacterium, so 
that the jump times become deterministic, and will occur on regular time instances 
separated by e 2 /A on the diffusive time scale; the randomness is then transferred 
to the jump size. The auxiliary process is thus defined as 

LA t7e 2 J-l 

X-t e = X c Q + ^ l+1 V n + (X t/e 2 - [X t/e 2 \)eAT[ Xoi/e2l+1 V [Xo - t/ ^ . 

n=0 s v ' 

ri 

{eq:sum_X_c} (4.6) 
The latter will be compared with a classical drifted random walk defined by the 
Markov chain 



{eq:Markov} { -"" + £ Aq +e X 2 \^n) v n, (47) 

v( 1+t ^ TRlv ')) v - 

to which we can associate a continuous process defined on a similar diffusive time 
scale by interpolation, 

/LAo*7e 2 J-l \ 



{eq:sum_Xi_c} (4.8) 
Now the computation of the diffusive limit of the sequence of processes 
(tt-^-X^' e ^j for e— >0 relies on four steps. 

Step (i) Show tightness of the coupled sequence [i <-> (X|' e ,S|) J . To this end, 

V / €>0 

decompose the leading terms in these processes into a "martingale" term 
and a drift term. 

Step (ii) Show, using a Gronwall argument, that the process describing the bac- 
terial position at jump times X^' e and the classical Markovian random 
walk have the same advection-diffusion limit. 
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Step (iii) Compute the advection-diffusion limit of the Markovian random walk S| 

(and hence of the process at jump times X^' e ). 
Step (iv) Use a random time change to transfer the result on bacteria position 

process X^' e itself. 

Step (i). Tightness of the coupled sequence t M> (X^' £ ,R^). 
Let us recall that according to Kolmogorov's criterion, a sequence of continuous 
processes (t h ^ w |) e>0 is tight in the space C([0,T],R d ) of continuous processes en- 
dowed with uniform convergence if there exists a constant C independent of e such 
that 



E\uf-ulf<C(t-sy 



(4.9) {eq:Kolmogorov} 



for any 0<s<t<T. In our case, checking (|4.9j) for the coupled process £i-> 



(Jf?' e ,Hj 6 ) is highly standard (see <^*). Remark that it is sufficient to check the 
criteria for each of the two separately. We recall the computation for the process 
fi— >S?' e , the case of £">— >X|' e being similar using (|4.5[) . First, one can check (|4.9I) 
when t — s< Aoe 2 . Indeed, in the latter case, only the rest r-i (see equation (|4.8[1 ) is 
involved and 

\4 



E 



4 <C^^e 4 <C(t-sf, 



where we have used that 



(4.10) {eq: est-eps-here} 



Assume now that t — s>Aoe 2 . The terms coming from r 2 can be easily controlled, 
since, again using (I4.10p . 

E|r 2 | 4 <Ce 4 <C(f-s) 2 . 

The terms of order e in the definition (|4.7[) (defining a martingale) can be controlled 
by expanding the sum and eliminating the terms with null average. It yields: 

4 



E 



E < 

ri=[Aos/e 2 J 



-iV B 



An 



<Ce 4 (LA f/e 2 j-LA s/e 2 j) 2 

<Ce i (i/e 2 -s/e 2 ) 2 + Ce i 
<C{t-sf<C{t-sf . 

In the last line, we use the fact that (i— s) 2 < 4T 2 < C. Finally the terms of order e 2 
in the definition (|4.7I) (defining a drift) can be controlled as 

4 



E 



|A i/e 2 J-l 

£ 2 v "-- 1 v •• i ./ 

™=L A oS/e 2 J 



\ 2 
A 



< 



Ce 8 (LA t7e 2 J-LAo5A 2 J) 4 



<Ce*(i/e 2 -- S /e 2 ) 4 + Ce* 
<C(i-sf, 
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where we again use the fact that (t - s) < AT 2 < C. We conclude that 

4 



E 



1 t 



<c(t-sy, 



and hence the sequence is tight. 

Step (ii). Comparison of the limits. 
We can now use a Gronwall argument to show that t^rX^ and i\- have the 
same limit (in distribution) when e— »0. First, since the coupled process is tight, we 
can use the Skorokhod embedding theorem, and consider a new probability space 
where an e-sequence converge almost surely (in uniform norm) towards a limiting 
process in > (X^'°,S?). Then, we compare the two processes using (|4.5j) and (|4.7|) . 
Incorporating the terms r\ (see equation (|4.6I) ) and r^. in the sum, we get 



< 



LA *7e 2 J 

E 

71 = 



\ 2 
A 



(A e (X^)-A e (E c Sj\+9 n+1 0(e 



so that taking expectation yields 

[\ i/e 2 \ 



E 



n=0 



Then, since the velocity of both processes is of Oe, and using the fact that A e is 
Lipschitz uniformly in e, we get 



(A e (Xy-A e (Z c n )) 
Finally we get 



<Ce- 2 



(n+l)e 2 Ao 



ne 2 /X 



x^ 



ds + e n+1 0(e). 



<c 



ds + 0(e), 



and applying a Gronwall argument yields 



E 



Xi 



<0(e 



Taking the limit e — >• shows that = S° almost surely. 

Step (iii). Compu tati on of the limit. 
Standard results (see ^3) show that the sequence of processes it- ► Sf satisfying (|4.7[) 
converges in distribution for the uniform norm towards the process solution to the 
SDE with diffusion matrix 



9 2 

7 n+ 



lV„<g>V„ 



2D 

V 



and drift vector field 



6(x):=A E 



'n+l 



V n A T (x)V r 



DA(x) 
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with D given by (|4.2[) . which yields the SDE (14.1[) . By Step (ii), the same holds true 
for t^X^' e . 

Step (iv). Random time change. 
It remains to show that t— > X^ and t— ^X|' e converge in distribution to the same 
limit. Let us introduce the random time change 

a c ,£(i) = e 2 T L c Aot - /e 2j +f 2 {T{ Xa t/e^\+l~ T {\ a i/^\) (V/e 2 ~ LV/e 2 j) (4.11) {eq:timechange_c} 

|A f/ e 2 J 

= £ 2 A^+ e 2 (Ar L c Aot7e2J+1 )(A t-A 2 -LAot"A 2 J), (4.12) 



which is exactly the linear interpolation of the jump times {T£) n>0 , and by con- 
struction is almost surely strictly increasing and continuous. Using (|4.5[) . we get 

2 l^oi/e 2 ] 3 [\ t/e 2 \ 

n—1 n—1 

so that 1 1 — y a c . e (t) converges in distribution on any time interval [0, T] for the uniform 
convergence topology towards the deterministic function (see e.g. ^J). By 

construction, 

X°^ ^=JTj' e , (4.13) {eq: alphatilde_c} 

so that we can apply the technical Lemma l4~2l to conclude. □ 



Remark 4.4. In one spatial dimension, one can describe the evolution of the 
distribution density p € (x,v,t) of the process with direct gradient sensing as fol- 
lows We introduce the distributions of left-moving and right-moving particles 
p±(x,t) — p c (x,v — ±1,£), and write, using ([33]): 



a c , a c A^(x,+1) A£(x,-1) 
+ e9xP+ = "P+ + £ P- 

q c q c A=(x,+1) AJ(x,-l) 
ftp! - ed x p c _ = p% - pi 



(4.14) {eq: cprocess_hydrold} 



If we assume X^(x,v) to be a linear function, 

\%(x,v) = A -eA(x)u, (4-15) {eq:lin-turn-lD} 

is easy to see that equation (|4.14[) is equivalent to 

d 2 n c + A d t n c = d x (e 2 d x n c ~ e 2 A(x)n c ) , (4.16) 

which, on diffusive time-scales i—te 2 , is equivalent to (|4.4I) in the limit when e— s>0. 
Note, however, that the assumption (|4.15[) requires |A(a;)| < Ao/e to retain positivity 
of the turning rate. For e small enough compared to A(x) (containing VS(i)) this 
requirement is always fulfilled. 
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4.2. Diffusion limit of the process with internal state 

In the same way, a standard probabilistic diffusion approximation argument can be 
used to derive the pathwise diffusive limit of the process with internal state (|2.5|) . 
We denote by 

the process with internal state on diffusive time scales. We will show that this 
process converges towards a solution of the same advection-diffusion SDE (|4.1[) as 
the process with direct gradient sensing (|2.10p for e— s>0, provided A(x) is chosen 
according to 

{eq:control_field} A (x) = b T lim - — ^ VS(x), (4-17) 

t->o AoT e +Id 

assuming that the limit exists. The parameters 6, r £ , and Ao in equation (I4.17P 
were introduced in (|3.2p ~ (12.41) as parameters of the process with internal state, and 
IdeR™ xn is the identity matrix. If we introduce the bacterial density as 

n(x,t)= / / p(x,v,y,t)M{dv)dy, (4.18) 

this implies that n satisfies (|4.4|) on diffusive time scales as e — » 0. 

The proof proceeds along the same lines as the proof of Proposition ^. 31 However, 
the derivation of an asymptotic expansion of the jump time differences becomes 
more involved. We need the following auxiliary definitions and lemmas : 

Definition 4.5. We denote by m:R— s-R nx " the function 

{eq:m(t)} m(t) := tr e - (id-e - *^ 1 ) r £ 2 , (4.19) 

whose derivative is given by 

m'(t) = T e (id -e" 4 ^ 1 ). 

Lemma 4.6. For all isR, we have 

\\m'(t)\\<t, 

as well as 

IK*)||<y. 

Proof. We have m(0) = m'(0) = and m"(t)=c~ tTc \ By Assumption 13.11 on r e , 
we have 

so that 

\\m'(t)\\ 



IK(*)II<1, 



m"{s)ds 



< I ds = t, 
o 
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as well as 



\\m(t)\\- 



m'(s)ds 



t 2 

< I sds = — . 

o 2 



□ 



Lemma 4.7. Let us denote (V 2 denotes the Hessian) as 

Cg = || V 2 ^!!^ , (4.20) {eq:nonlin_signal} 

as well as the error term due to non-linearities of the problem as (recall that by 
assumption 5 <l/k) 

r (e) = c x e kS +c F e 1+S + c s e 2 <0(e), (4.21) {eq:nl_err} 

where the constants are defined in resp. equations (12. 4j) . (|3.2p . and (|4.20j) . The 

ference between two jump times, 



AT„ +1 :— T n+ i — T n 

satisfies the equation 

O n+1 = X AT n+ i~b T m'{AT n+ i)Z Tn 



- c 



(4.22) {eq:DTl} 

b T m(AT n+1 ) VS(X Tn )V n + (0 3 n+1 + 9 n+1 ) O (r(e)) . 



Proof. Throughout the proof, we will use the estimate |AT n+ i| < C9 n +i : with- 
out always mentioning explicitly. Duhamel's integration of the ODE in (|2.5p on 
[T n ,T n+ i] yields 

Zt=e -(t-T n )r^ ZTn + e C e -(t-t')^- l V5(JCtOVt/d*' 

nt 

+ e~( t ~ t ') Tl ~ 1 (r c _1 Z t / — F t (Y t i , S (X t ')) dt! . (4.23) {eq:Duhamel_lin_2} 

Since X t , = X Tn +e(t'-T n )V n , with t' <T n+1 , we get, 

| VS(X t , ) - VS(X Tn ) | = 9 n+1 c s O(e) , 

using cs as defined in (|4.20[) and the estimate t' — T n < AT„ + i < C8 n +±. Moreover, 
Assumption 13.11 combined with Assumption 13. 4| ensures that 



\T £ - 1 Z t ,-F e (Y t ,,S(X t ,)\=e 1 - s c F 0(\Z t ,\ 2d ) = c F 0(e 1+s ). 
Inserting these two estimates in (|4.23p . yields 

Z t = e- {t - T ^ 1 Z Tn +e J e-( t - t 'h~ 1 dt'WS(X Tn )V Tn + (.0l +1 + e n+1 )O(c F e 1+s + c s e 2 ) 

= e -(t-T n )r^ ZTn+e ^ Id _ e Ht-T n )r^y eVS{XT ^ 

(4.24) {eq:key_Duhamel} 
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Integrating again, we get 

J " +1 Z t dt= (id - e~ (T " +1 " T " )Te " 1 ) r e (Z Tn - er e VS(Xr„ ) V T „ ) 

+ e(T n+1 - T n )T f VS{X Tn )VT n + (0f l+1 + 6l +1 )0(c F e 1+s + c s e 2 ) 
= m\AT n+1 )Z Tn +em(AT n+1 )VS(X T JVT n + (6l +1 +6 2 n+1 )0(c F e 1+s + cse 2 ), 

with m(i) and m!(t) denned as in Definition 14.51 Finally, (|2.4[) gives 



T, 



+ i 



n+1 = A o AT„ +1 -fc J / Z t a!i + #„ +1 c A 0(e M ), 

and the result (|4.22|) follows. □ 

The following estimate of jump times then follows : 
Lemma 4.8. The jump time variations can be written in the following form : 
{eq: DTestim} AT n+1 = AT° +1 + eAT,J +1 + {9 6 n+1 + 6 n+1 )0 (r(e) + e 2 ) , (4.25) 



where 



{eq:estimDTO} AT„° +1 = + ^m'(AT„° +1 )Z T „ , (4.26) 



{eq:estimDTl} AT* +1 = -— 6 T m(AT° +1 )VS(X T JV„, (4.27) 

or, /or short, 

{eq:estimDTO_bis} AT% +1 = +6 n+l O(e S ) (4.28) 

{eq:estimDTl_bis} AT,\ +1 = -U T m ( -2±i ) V5(X T „)V„ + 6» 2 + 1 C(e' 5 ), (4.29) 
where the term of order e s in (|4.28|) is independent of V n . 

Proof. Let us first introduce the functions 

(y Q (t) = \ t-b T m'(t)Z Tn 
\y 1 (t) = -b T m(t)VS(X T jV n , 

Using these functions, we can rewrite the result (|4.22j) of Lemma T4.71 as 

{eq: DTlbis} n+1 = ¥„( AT n+1 ) + e*i(AT n+1 ) + (0, 3 1+1 + 0„+i) O (r(e)) . (4.30) 

Also, the definitions (|4T26]) - (|4T271) of (AT° +1 , AT^ +1 ) can be written as 

tt o (AT° +1 ) = n+1 

{eq:psi_DT} <J ^(AT° +1 ) (4.31) 

- +i n(AT° +1 )- 
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Taking the first line, add e times the second line, and rearrange, one has 

{e q :why_this} 9 n+1 = * (AT° +1 ) + e%(AT? l+1 )AT,\ +1 +e* 1 (AT° +1 ). (4.32) 

We then consider the expansion 

AT n+1 = AT„° +1 + eAT,i +1 + R n+1 . 

The proof of the Lemma then reduces to showing that R n +\ is small. Anticipating 
the result, the idea is to Taylor expand $0(^+1) and $i(AT, 1+ i), use these 
Taylor expansions in (|4.3ip . and compare to (|4.32l) . We start by using a double 
Taylor expansion at order 1 for VE" o(AT„ + i) : 

* (AT n+ i)=*o(AT° +1 +i?„ +1 ) + eAT I i +1 ^(Ar° +1 + i?„ +1 ) 

+ ^ +l)2 <(Ar»° +1 +Rn+i +73eA^ +1 ), 
=*o(ATO +1 ) + J R n+1 ^(AT° +1 +7 1J R n+ i) + eATi +1 ^(AT° +1 ) 

+ eAT„ 1 +1 i?„+i*d'(AT^ +1 + 72 i? Il+1 ) + <(AT° +1 + + 73 6AT„ 1 +1 ), 

and a simple Taylor expansion at order for ^i(AT n _|_i): 

* 1 (AT n+1 ) = * 1 (Ar° + + (eAr n 1 +1 +E„ +1 )^(ATO +1 +74(eAT n V + iin + i)), 

where 7$ € [0, — 1 . . .4 are used to apply the intermediate value theorem to the 
rest of the Taylor expansions. Inserting the Taylor expansions in the identity (|4.30j) . 
and comparing to (|4.32l) simplifies to 

0=%(AT° +1+11 R n+1 )R n+1 

( AT" 1 ) 2 

+ eATi +1 R n+1 V>(AT° +1 +72^+1) + ^— f ± ^-*o(AT° +1 +i?„ +1 + 73 eAr n 1 +1 ) 

+ (e 2 AT n x +1 + ei? Il+1 )*i ( AT° +1 + 7 4 (eAT^ +1 + i?„+i )) + + n+1 ) O (r(e)) , 

(4.33) {eq:DTlter} 

We now need some estimates on the terms involved in (|4.33|) to control R n +i- 
First of all, we need some crude bounds on jump times. Since rates are bounded 
above and below by A m i n and A max , we have 

\AT n+ i I < C8 n+ i. 

Then, *&' (i) = Ao — \ob T e~ tTe 1 Zr n , so that, using Assumption I3.4[ there is a con- 
stant C such that 

1 3 
-A < sup %{t)<-\ 0) 

i «e[o,+oo[ L 

for e < I/O. As a consequence, AT° +1 is well-defined by the implicit definition (|4.26l) 
and one has obviously 

AT n+1 <C8 n+1 . 
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In the same way 

\ATi +1 \<Cm(AT° +1 )<C6 2 n+v 

and thus by definition: 

< C(6 n+ i + e9 2 +1 ). 

Then we need some bounds on the derivatives of V^o and ^i- Using Lemma I3T41 we 
get 



6 T rr 1 e- tr ' 'Z T 



<C. 



\<(t)\ = X 
Eventually, for any t > 0, 

Wi{t)\<Cm'{t)<Ct 2 . 
Inserting all the previous estimates in (|4.33l) to control R n +i, we get 

\Rn+i\ < Ce(0i +1 + 9 2 n+1 )+(d 3 n+1 + e n+1 )O(v(e)) . 
Inserting a second time in (|4.33[) to control Rn+i, it yields: 

\R n+ i\<Ce 2 (0 6 n+1 +el +1 ) + (ef l+1 +e n+1 )o(T(e)), 

and the result follows. □ 

Remark 4.9. It is crucial to note that the term of order e s in (|4.26|) - (|4.28|) is 
independent of V n . This will ensure that the latter contributes to martingale terms 
in sums of the form ^2 n ^T n V n . 

Finally, the drift of the limiting diffusion (see Proposition below) will be com- 
puted using the following lemma : 

Lemma 4.10. Let 6 be an exponential random variable of mean 1. Then: 

<•*-■-> E K^)Hra (434) 

Proof. It is a straightforward computation: 

n + OO 



E ( TO (x")) = / rn{t/X)e' t dt 



= IL / te-'dt-r 2 / e-'dt + r 2 I ' e-^+^'dt 



Ao 

:^-T e 2 +T e 2 (ld+T £ -7A 



-1 



1 2 



Id+A T e \X 
A Id+A T, 



— (Id+A r e ) (Id-A r e ) + A t £ 
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□ 

Using these lemmas, we can prove the following result on the diffusive limit of 
the process with internal state : 

Proposition 4.11. Consider the process with internal state define in Section \2.1[ 

and the assumptions of Section \3.3l Assume lim^o t r exists, and that 

A T e +id 

5>l/k. Then, the process t<-^X^ = X t / e 2 converges in distribution (for uniform 
convergence topology) towards i>—>X® solution to the SDE (|4.ip with drift coef- 
ficient (|4.17p . As a consequence, the density nf of the process converges weakly 
towards = solution of {4-4ty with drift coefficient (|4.17[) . 

Proof. Since the proof is very similar to the proof of Proposition ^. 3[ we will skip 
similar arguments and refer to the latter whenever possible. As in the proof of 
Proposition 14.31 we perform a random time change on the process of positions of 
the bacterium, so that a jump (tumble phase) become deterministic, and will occur 
on regular small time interval of length e 2 /Ao. The latter is thus defined as: 

|A t/e 2 J-l 

Xf = X + V eAT n+1 V n + (X i/e 2 -lXoi/e 2 \)eAT lXoi/e 2 i+1 V^ oi/e 2 l . 

n=0 " - ' 

ri 

(4.35) {eq:sum_X} 

To compute the diffusion limit X| , we will compare it to a classical drifted random 
walk defined by the Markov chain 

„ a (4.36) |eq:Markov int) 

To the latter, we can associate a continuous process defined on diffusive time scales 
by interpolation 

/l\ n i/e 2 \-l \ 

£f = Jf + ^ s n+i-S n +(X i/e 2 ~[X i/e 2 \)(Eix o i/^\+i-^[\ot/^\)- 
I n=0 / v • ' 

(4.37) {eq:sum_Xi} 

Now the computation of the diffusive limit of the sequence of processes i >-> X$ relies 
again on four steps. 

Step (i) . Tightness of the coupled sequence t M> (X| , H|) . 
The proof of tightness is strictly similar to the proof in Proposition 231 so we leave 
it to the reader. 

Step (ii). Comparison of the limits. 
We can now use a Gronwall argument to show that t\->Xi and f^S^ have the 
same limit (in distribution) when e— >0. First, since the coupled process is tight, we 
can consider a new probability space where a e-sequence converges almost surely 
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(in uniform norm) towards the limiting process ti— > (Jr5,E|). Then, we compare the 
two processes using fl4.25[) - (|4.35|) and (|4.36[) - (|4.37p . For this purpose, we incorporate 
the rest terms r\ and in the sum in (|4.35|) and (|4.37[) . and we get the estimate 



Xf -El < 

|A t/e 2 J 



n— ^ ^ 



Ao Aq 



+ J2 ^ 2 V„ A^ +1 --m(^±i)V5(S„)V„ + E (^ +1 +^+i)0(e 3 + er(e)). 

„= v A A y n=0 



Denoting 

and using the estimates in Lemma 14.81 we get 



M n = V n b T m'(AT? l+1 )Z Tn , 



XI -Sf 



< 



LA t7e 2 J 

E eM " 



|a *7« 2 J 

E 



n=0 



2 °n+l 



A 







Ao 



LA *7e 2 J 



|VS(X r *)-VS(S n )| + E (0 n+1 +e n+1 )O(e 3 + e 2 + s + er(e)). 



71=0 



Taking expectation and using Jensen's inequality, as well as the fact that S is 
Lipschitz, then gives 

1/2 











LAo*7e 2 J 




E 


xi -n 


< 


E 


E eM n 










V 


n=0 





+ Ce 2 E IE|X T c-S n |+0(e + e 5 + e- 1 r(e)). 

n=0 

Using \M n \ <C9 n+ ie s (which follows from Assumption 13.41) and the independence 
of the random variables (Vn) n >o, we get 

2 

<C[X t/e 2 \ xe 2 xe 2S = 0{e 2S ). 



E 



LA *7^J 

E eM - 

n=0 



Controlling the derivative of both processes on the time interval [ne 2 / Xo,(n- 
l)e 2 /A ] yields 











Xt% — S„ 


<A e- 2 / 


x g — S| 


ds + n+1 O(e) 




J ne 2 / Ao 







In the end, we get 

E Xf-Sf 



<C / E 

o 



Xg — 3| 



ds + C(e + e <5 + £- 1 r(e)) 
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and applying a Gronwall argument yields (expliciting r(e)): 



{eq: f inal_bound} E 



■■0((l + c s )e+(l + c F )e d + c x e k5 - 1 ). (4.38) 



Taking the limit e — > shows that X^ = almost surely. 

Step (iii). Computation of the limit. 
Standard results (see^) show that the sequence of processes i >->• S| satisfying (|4.36|) 
converges in distribution for the uniform norm towards the process that is solution 
of the SDE with diffusion matrix : 

A., .MA, 
and drift vector field, using (|4.34|) , see Lemma 14.101 : 



6(a;) = limAoE ^b 1 m(^WS(x)V r , 



= limD(6 T — ^— VS(x))/X , 
m*DA(x)/\o, 



with A(x) defined in equation 14. 171 yielding the SDE (|4.1j) . By Step (ii), the same 
holds true for ii-^X}. 

Step (iv). Random time change. 
It remains to show that i— >X^' e and i—^X^ e converge to the same limit. We 
introduce the random time change 

a e (f) = e 2 T LAof/e 2j +e 2 (Ty Xo t/^\ + l-Ty Xo t/ e 2^) (A i/e 2 - \\oi/e 2 \) (4.39) {eq:timechange} 

LA t7<s 2 J 

- J2 € 2 AT n + e 2 (AT lXo t /e 2 l+1 )(X t/e 2 ~[X t/e 2 \) : (4.40) 

which is the linear interpolation of the jump times (T n ) n>0 , and by construction 
is almost surely strictly increasing and continuous. Lemma 14.81 implies that AT n — 
2 ^T + (Qf l +i+Qn+i)0(e + e 1+s + r{e))). Therefore, we get 

«.©=f E + f +£l(£)) E (#+«Wi)> 

n— 1 n— 1 

so that <i— > a C!£ (i) converges in distribution on any time interval [0,T] for the uniform 
convergence topology towards the deterministic function ti— >i (see, e.g., ES). By 
construction, 

X k{i)= X h ( 4 - 41 ) {eq:alphatilde} 

so that we can apply the technical Lemma |4~21 to conclude. □ 
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Remark 4.12. In one spatial dimension and with internal dynamics given by (|2.8p . 
Erban and Othmer argued on formal arguments (a moment expansion) that the 
evolution of the density n(x,t) on hyperbolic time-scales satisfies 

d^n + X^dtn — dx ( e 2 d x n— ^—S'(x)n], (4.42) {eq:process_hydro} 

V 1 + A t £ / 

in the limit when e->0^. It is easy to see that equation (|4.42|) is equivalent to (|4.1j) 
in this limit when choosing A(x) in (|4.ip according to (|4.17[) . 

Remark 4.13. Another relation between the process with internal dynamics and 
the process with direct gradient sensing is obtained when, for fixed e, the bacteria 
learn the value of the chemoattractants concentrations infinitely fast. To investigate 
this limit, we may construct a sequence (T e fe )^° =1 , such that lim^oo ||r e fc || _1 = oo, 
and show that, for k— >oo, the process with internal state converges to a process 
with direct gradient sensing, with jump rate given by A°(x,v) — X e (VS(x)v), i.e. 
the turning rate for the process with internal dynamics, evaluated for the deviation 
z = VS(x)v. 

5. Discretization of velocity-jump processes 

To observe the asymptotic behaviour in numerical simulations, the time discretiza- 
tion should approach the same diffusion limit as the time-continuous problem. This 
requires, in particular, that the discretized jump times are computed sufficiently 
accurately. In this section, we propose a time discretization of (|2.5|l and give a con- 
sistency result. The process with direct gradient sensing (|2.10l) can be discretized 
similarly. For ease of exposition, we consider the scalar equation (|2.8p for the internal 
state; generalization to nonlinear systems of equations is suggested in Remark [37TJ 

5.1. Linear turning rate 

First, assume a linearized jump rate 

{eq: jump-lin} \(z) = \q — b T z. (5-1) 

When the chemoattractant profile is linear, S(x) = So + (V5o) ■ x, the ODE (|2.8[) 
can be explicitly solved for Z t . For i€ \T n ,T n+ i), i.e. between two jumps, we then 
have 

( X t = X Tn + eV t (t-T n ) 

{e:discr-lin} < _ t / T -i\ v (5.2) 

|z t = e-(*- T "^ Z Tn +eT e (ld-e-^- T ^ )(V5 )V„. 

Similarly to the computation in the proof of Lemma 14. 7[ we get, using AT n+ i :— 



\(Z t )dt = \ Q (AT n+1 )-b T (ld-e- AT ^ T <~ jr e Z Tn 
{e: jump-lin} - eb T (AT„ +1 T e - (Id -q-^+it.- 1 (y Sq ) V„. (5.3) 
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The jump time T n+ i can be calculated exactly from (|5.3[) by solving 

" +1 \(Z t )dt = e n+ll (5.4) 

using a Newton iteration. Hence, no time discretization error is made. 

When S(x) is a nonlinear function of x, exact integration of Yt is no longer 
possible. We therefore define a numerical solution (Xf* ,V t St ,Y t St ) as follows. Between 
jumps, we discretize the simulation in steps of size St and denote by {Xff kl Zff k ) 
the solution at tn^ — Tf^ + kSt. The numerical solution for t£ [t n: k,t n ,k+i] is given 

by 

(x^=x^ k + e v n {t-t n , k ) 

{Zp =exp(-(i-^ fc )T e - 1 )Z5 + er e (ld-exp(-(t-t n , fc )T e - 1 ))V5(X I 5 l * fc )V„. 

(5.5) {e :discr-nonlin} 
We denote by K > the integer such that the simulated jump time T^ +1 £ 

rpSt 

[t n ,K,t n ,K+i]- To find T^ +1 , we first approximate the integral f T s t +1 \(Zf t )dt us- 
ing 

f " +1 A(Zf )dt= £ f "' fc+1 A(Zf )dt + f " +1 A(Zf )dt, (5.6) 

JT&* k=0 Jt n , k Jt ntK 

and then compute: 
pi 



n k+1 X{Zf t )At = X 5t-b T [ld-e- 5tT ' l )nZfi n 

-eb T (Str e - (Id-e-^" 1 )^ 2 ) VS(X% k )V n . (5.7) 
The jump time T^ t +1 can then be computed as the solution of 

/*»,*+! 

/ A(Zf )dt = 6> n+ i- 2J / ^(Zf)dt, (5.8) {e:lin-Newton} 



again using a Newton procedure. The use of (|5.5[) as a discretization scheme has 
a crucial consequence: the analysis performed in Section [3] for the continuous time 
case still holds in an exact fashion in the time discretized case. This is mainly due 
to the fact that the key estimate (14.24|) still holds. It is left to the reader to check 
that all the proofs performed in Section [4] are not affected by this discretization. 

Remark 5.1. The discretization scheme (|5.5[) can be generalized to a non-linear 
evolution equation for the internal state. It is however necessary to check that a 
time discrete version of Assumption 13.41 still holds in this case, and also that (|4.24[) 
remains true. 

We also give the following consistency property, which is of independent interest : 
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Lemma 5.2. Assume T^ = T n , X Tn =Xft n , V n = V^ and Y^=Y Tn , then, using 
the time discretization (|5.5p . we have 

.2 r,/)2 



and consequently, 



\T n +i -Tn + ± \ < Ce z 5td^ +1 , 



Xt„ +1 — X!^s t <Ce 3 St9 2 l+1 . 



(5.9) 

(5.10) 
(5.11) 



Proof. Using the assumption T n = and Xt„ = Xjf , we get straightforwardly 



Xt=X° l , yte[T n ,T n+1 ATX +1 ] 



so that 



vs(x t )-vs(x s n t 



< CSte. 



Consider Duhamel's formula (|4.23|l . and write a similar expression for the time 
discretized case for 1 6 [T n ,Tn +1 ] : 

rt 



{eq : Duhamel_num} 



{eq:dZ} 



)/**J 



V n dt'. (5.12) 



Using the assumption = Zt„ , the difference between the two Duhamel's integrals 
then yields 



\Z t -Z d t t \<Ce 2 Std n+ i, 

ioTte[T n ,T n+1 AT^ +1 }. 

We now proceed to find a bound on \T n+ i —T^ +1 \. We start by writing 



(5.13) 



A(Z t )dt = 



T„ 



A(Zf)di 



and using the assumption T n =T^ t , we obtain 



L 



b T (Z t -Z d t t )dt = 



A(Zf)di 



T„ +1 AT*^ 1 



A(Z t )dt 



T„ +1 AT*^ 1 



J t% ) 



= (T^ +1 ~T n+1 AT^ +1 ) A(Zfi) + (T„ +1 -T n+1 AT„« J A(i?| 

for some values of *J £ [T n+1 AT^ +1 ,T^ +1 ] and *| £ [T n+1 AT^ 1( Tf]. Since A min < 
AO?) < A max , we get 



rjiSt rji I x" v 

_J n+l St> 



T„ + iAT; 



Zt-ZfHdt 



<C9 2 n+1 e 2 St. 



Finally, note that 



Xt-Xp \<e(t-T n ) 
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for t>T n to get the estimate on 



-X. 



at 



This concludes the proof. 



□ 



Note that the dependence in e 3 5t shows, at least at the consistency level, that 
the error due to discretization remains of a higher order in e than the difference 
between the process with internal state and its advection-diffusion limit. 



5.2. Nonlinear turning rate 

We now consider the nonlinear turning rate as defined in (|2.9j) . We again discretize 
in time to obtain the time-discrete solution (|5 . 5|) . The jump time T^ +1 is now 
computed by linearizing (|2.9|) in each time step, 



A* (Zf* A) = A(Z* )+ dA(Z " tfc) ^ 



dz 



r/Ot ryo 



(5.14) {eq:nonlinrate_disc} 



and approximating the integral 



rpSt K—l 



\ St (Z*\Z s n %)dt + / \ St {Zt\Z s *)dt 



7 n+l, 



(5.15) 



The jump time T^ +1 can then again be computed using a Newton procedure as 
in (|5.8p . Here again, the analysis performed in Section [4] for the continuous time 
case, still holds in an exact fashion in the time discretized case. Also, a consistency 
property similar to Lemma 15.21 holds in the case of non-linear turning rate. 



Lemma 5.3. Consider a non-linear turning rate (I2.9[) discretized by (|5.14[) and 

| . Assume T r f = T n , X Tn = Xf? , V„ = Vft and = Y Tn , then 



T n +i — T^+i | — Ce z (St9f l+1 +St /: 8 n+ i) 1 



and consequently, 



Si 

rp5t 



<Ce 3 {St9l +1 +St 2 6 n+1 ). 



(5.16) 



(5.17) 



Proof. Similarly to the linear case we can show that 



Irp rriSt I ^ S~1 

\ 1 n+l ~ J- n+l -i <- / 



t„+iat; 



+l AJ n+l 



\{Z t )-\^\Zt\Z 6t t ^ t 



dt. 



We then estimate the righthand side of (|5.18j) . This is done in two steps 



(5.18) {eq:diffl} 



Step (i) Estimate difference between turning rate and locally linearized turning 
rate for the discretized deviations A(Zf ) - A 5t (Zf ,Z** fe ) . 

Step (ii) Estimate difference between turning rates for the discretized and time- 
continuous deviations \\(Z t ) — X(Z^ t )\. 
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Step(i): \{Z?)-\St{Z s t \Z% k ) . 
A Taylor expansion for t£ [t n ,k^n,k+i\ yields 

| A(Zf ) - \ St (Zf ,Z% k )\<C {Z? - Zi%) 2 . 

Then, remarking that, for tE [t n ,k,t n .k+i], we have 

^ = e-( t -*-.^«- 1 Z« fc + cr e (M-e(-( t -*-.^- 1 ))v5(^) ) 

we get 



I iySt iySt I . 

\ Zj t ~ Zj nM~ 



Id — exp I — 



t t r . 



{Z 5 * -r e eV5(X" )) 



= |m'(t-i n , fc )(ep S V^(^ fe )~r e - 1 z5)| 
< CJie. 

In the last line of the above, we have used the estimate 

\Zf \=0(e 5 ) 



as in the non discretized case, which follows from Duhamel's formula (|5.12|) and the 
assumption I —0{e s ). Finally, we obtain 



|A(zf)-A 5 *(Zf,^ fe )|di<Cft 3 



and thus: 



\{Z t )-X>*[Z$\Z>* tBi0 _ t 

n > L St J 



dt<C<5i 2 e 2 6»„+i 



Step (ii): \\{Z t )-\{Zf% 

Using Duhamel's formula as in the proof of Lemma 15.21 and using (|5.13[) above 
yields 

I \{Z t ) - A(Zf ) I < C \Z t - Zl l I < Cfe 2 0„ +1 . 



Finally we get, using (|5.18JI . 

\T n+1 -T^ +l \ < Ce 2 (6t6l +1 +5t 2 9 n+1 ). 
This concludes the proof. 



□ 



6. Conclusions and discussion 

This paper proved rigorously, using probabilistic arguments, the convergence of a 
velocity-jump process for chemotaxis of bacteria with internal state to an advection- 
diffusion limit, and showed that the same advection-diffusion limit can be obtained 
as the diffusion limit of a simpler, coarse process with direct gradient sensing. From 
a numerical point of view, we introduced a time-discretized model that retains these 
diffusion limits. 
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A direct simulation using stochastic particles presents a large statistical variance, 
even in the diffusive asymptotic regime, where the limiting behavior of the bacterial 
density is known explicitly. Consequently, it is difficult to assess accurately how 
the solutions of the fine-scale model differ from their advection-diffusion limit in 
inter medi ate regimes. Therefore, numerical results are deferred to the companion 
paper ^H, where we will use the results of this paper to construct a hybrid scheme 
for the simulation of the fine-scale process with "asymptotic" variance reduction. 
The hybrid scheme simulates both the process with internal state and the simpler 
process with direct gradient sensing using the same random numbers, and exploits 
a pathwise coupling between the two trajectories. The variance reduction that is 
"asymptotic" in the sense that the variance vanishes in the diffusion limit. 
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